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1 Introduction 

One of the questions which ought to be addressed before performing a large 
scale computer simulation is "What are suitable weight factors for the problem 
at hand?" It has been expert wisdom for quite a while that Monte Carlo (MC) 
simulations with a-priori unknown weight factors are feasible and deserve to 
be considered 1 . With focus on narrow classes of applications, this idea was 
occasionally re-discovered, for instance 2 . But it needed the work of Ref. 3 ' 4 to 
become a more widely accepted idea in Statistical Physics as well as in Lattice 
Gauge Theory. 

The next section introduces systems with supercritical slowing down. Sec- 
tion 3 reviews the multicanonical (MUCA) approach and briefly sketches re- 
lated methods. Applications to first-order phase transitions are summarized 
in section 4 and those to systems with conflicting constraints in section 5, 
before some conclusions are drawn in section 6. Some interesting topics, like 
the dynamical-parameter method in U(l) gauge theory 5 , had to be omitted 
because of space limitations. For a MUCA calculation of constraint effective 
potentials see Neuhaus in these proceedings. 



2 Supercritical Slowing Down 

The terminology "supercritical slowing down" is used to characterize canonical 
MC simulations which slow down exponentially fast with increasing system 
size. It is useful to distinguish static and dynamic reasons for the slowing 

"Published in the Proceedings of the International Conference on Multiscale Phenomena 
and Their Simulations (Bielefeld, October 1996), edited by F. Karsch, B. Monien and H. 
Satz (World Scientific, 1997). 

^E-mail: berg@hep.fsu.edu; WWW: http://www.hep.fsu.edu/" berg. 



1 



down. In the former case desired configurations have an exponentially small 
weight, whereas in the latter case their weight may still be large, but the 
dynamic process of reaching them deteriorates. 

Static Examples: 

• Magnetic field driven first-order phase transitions: Configurations with 
zero (or small) magnetic fields are exponentially suppressed at low tem- 
peratures and they exhibit domain walls. 

• Temperature driven first-order phase transitions: Configurations with 
domain walls are exponentially suppressed. 

Dynamic Examples: 

• Low temperature transitions between magnetic states (for instance the 
up-down states of the Ising model below the Curie temperature, ...). 

• Transition between low temperature states in systems with conflicting 
constraints: Spin glasses, the traveling salesman problem, ... 

3 Multicanonical and Related Methods 

These methods try to overcome supercritical slowing down by sampling with 
unconventional weights. 

Canonical MC simulations perform importance sampling with respect to 
the Boltzmann weight 

w B = e-P E 

Using reweighting, expectation values in a vicinity (— > for V — > oo) of the 
temperature T = 1/(3 are obtained. 

Multicanonical refers to simulations which obtain expectation values for a 
temperature range, which stays finite in the limit V — > oo. 

Similarly, multimagnetical 6 refers to simulations which give results for a 

certain range of the magnetic field, etc The cluster variant 7 is called 

multibondic. 

3.1 How to get MUCA results? 

The aim is to sample a broad energy density, like the uniform 
P(E) = const, for E min < E < £ max 
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where e m ; n = E m i n /V and e max = E max /V may be kept constant. The uniform 
density is obtained by sampling, for E min < E < E max , with the weight factor 

w(E) = l/n(E) 

where n(E) is the spectral density. As n(E) is a-priori unknown, some pre- 
liminary or iterative estimate of n(E) has to be part of the approach. For 
attempts to optimize the weight factors see 8 . 

3.2 Weight factors and temperature 
Let us re- write the weight factor as 

w(E)=e- s ^ = e -P(E)E +a{ E) 

where S(E) is the microcanonical entropy. Then, the temperature follows from 

T(E) = l/[3(E) and (3{E) = ^1 
Further, the function a(E) is determined up to an additive constant. 
Example with discrete energy: 

13(E) = [S(E + e) - S(E)]/e (e smallest energy step) 
The identity S(E) = f3{E)E - a(E) implies 

S(E) - S(E - e) = j3(E)E - (3{E - e)(E - e) - a(E) + a(E - e) 

Inserting ef3(E - e) = S(E) - S(E - e) yields the recursion 3 ' 4 

a(E - e) = a(E) + [f3(E - e) - f3(E)]E 

where a(E max ) = is a convenient choice of the integration constant. 

3.3 Recursive weight factor estimates 

For spin systems with first-order phase transitions the finite size scaling (FSS) 
behavior is relatively well-known. Provided the steps between system sizes 
are not too large, it is then possible to get working estimates of the MUCA 
weights by means of FSS extrapolation from the already simulated smaller 
systems 6 . Another method which works well for these systems is patching 
of overlapping, constraint 9 MC simulations. While these approaches seem to 
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work well for static slowing down, they fail for the dynamic slowing down of 
spin glass simulations. Then it is recommended to employ a more sophisticated 
recursion 10 , as outlined now. 

It is advisable to start a recursion for MUCA weights in the disordered 
region (for which reliable canonical calculations can be performed). For in- 
stance, 

$>{E) = /3°(£ max ) e.g. /?°(i? max ) = 

Then 

f3 n+1 (E) = e- 1 ln[HZ(E + e)/H$(E)] (1) 

where H™(E), (x = 0,(3) contains combined information from the runs with 
P°(E),...,p n (E): 

n 

H:{E) = Y,9 k {E)H k x {E) 

k=0 

Here Hq (E) is the unnormalized histogram obtained from the (short) simula- 
tion with (3 k (E). Further, 

H k (E)=H k (E) e~P k ^ e 

and the factors g k (E) weight the runs suitably. For instance, 

g k (E) = const. min[H k (E + e), H k (E)} 

where the constant follows from the normalization ^2 k g k = 1. For energy 
regions in which the statistics is (still) insufficient, equation (1) may be sup- 
plemented by 

f3 n+1 (E) = (3 n+1 (E + e). 

An interesting recent idea 11 is to utilize transition frequencies instead of simple 
histograms. Finally, connections with adaption and linear response theory have 
been explored 12 . 

3.4 Slowing down 
Our typical situation is 

Em&x -^min ^ U 

The MUCA optimum for a flat energy distribution is given by a random walk 
in the energy. This implies a CPU time increase 

~ V 2 
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to keep the number of E max — > E m i n — » _E max transitions constant. The 
recursion (1) needs an additional ~ V 0,5 (optimum) attempts to cover the 
entire range. It follows: 

slowing down ~ V 2 ' 5 or worse. 

Recursion alternative (patching of overlapping constraint 9 MC simulations) : 

number of (fixed size) patches ~ V . 

When results can be obtained by keeping the number of updates per spin 
(sweeps) in each patch constant, another CPU factor <~ V follows. In this case 
we can get: 

optimal performance ~ V 2 . 

This is still the optimal slowing down when the MUCA parameters can be 
obtained via FSS extrapolations. 

3.5 Related methods 

Combining MUCA with multigrid methods has been explored in Ref. 13 and 
the cluster version 7 has already been mentioned. For molecular dynamics, 
Langevin and hybrid MC variants see Ref. 14 . 

The method of expanded ensembles 15 proposes to enlarge the configuration 
space by introducing new dynamical variables. In simulated tempering 16 it is 
the temperature. A discrete set of weight factors is introduced 

w k =e-^ E+a ", k = l,...,n, ft < fa < - < Pn-l < Pn 

The transitions 

{Pk,ak) — » (/3k-i,afc_i), (P k +i,a k+ i) 
are now added to the usual E — > E' transitions. 

Remarks: 

• The method works for dynamic but not for static supercritical slowing 
down, because each member of the discrete set of weight factors sam- 
ples still a Boltzmann distribution (e.g. it is not suitable for calculating 
interfacial tensions) . 

• In the context of massively parallel computer architectures (and beyond) 
the variant of parallel tempering is particularly promising 17 . For addi- 
tional information see the article by Marinari in these proceedings. 
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Random Cost 18 sacrifices the exact relationship with the canonical ensem- 
ble in favor of having a-priori defined transition probabilities. Assume, we can 
choose from a discrete set of updates, such that an update implies one of the 
following energy changes: 

AE+, (i = l,...,n+), A£°, (i = l,...,n°) or A££, (fc = l,...,n") 

where A £7+ > 0, A_E° = and A_E,7 < 0. It is then easy to define update 
probabilities pf, p" and p^, J2iPt + Y,jP°j + J2kPk = sucn tnat 

^p+A^ + = -E^ A ^- ( 2 ) 

i k 

holds (but at the extrema). The algorithm performs a random walk in the 
energy and, hence, samples a broad energy distribution. It is conjectured to 
be of advantage in optimization problems, where one is mainly interested in 
minima and less in the canonical ensemble. 

4 First-Order Phase Transitions 

MUCA simulations are best established for investigations of first-order phase 
transitions. The range of applications goes from MUCA studies of mathemati- 
cally ambitious topics to chemistry oriented ones. A MUCA investigation 19 of 
the Borgs-Kotecky 20 FSS theory shows that very strong phase transitions or 
very large lattices are needed to observe the asymptotic behavior. To give two 
examples from the chemistry side, an investigation of the coexistence curve of 
the Lennard- Jones fluid was performed in Ref. 21 and the liquid- vapor asymme- 
try in pure fluids was studied in Ref. 22 . 

Most work has focused on calculations of interfacial tensions and an over- 
view is given in the forthcoming. 

4.1 2d Potts Models 

A pioneering MUCA study was performed for the 2d ten-state Potts model 4 . 
Using the histogram method 23 , the value 2/ s = 0.0978(8) was found through 
FSS study of the equation 

2f[ = -ylnPf a where Pf n = mm{P L (E) \ E™*' 1 <E< £™ ax > 2 } . (3) 
E 

Here Pl(E) is the normalized energy density at the temperature defined by 

pmax.l = p max^ when ^ MUCA estimate 

was published, numerical 
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estimates of Potts model interfacial tensions disagreed up to one order of mag- 
nitude. But, shortly after the exact value was found to be 2/ s = 0.094701..., 
see 24 ' 25 and references therein. Certainly, this helped towards the break- 
through of MUCA methods. Once, the exact result was known, the remaining, 
small discrepancy could be eliminated by improving the finite volume estima- 
tors 26 . 

For these simulations the MUCA slowing down is around ~ V 2 ' 3 , i.e. 
reasonably close to the optimal performance. For related investigations of the 
seven-state 2d Potts model, see Ref. 27,28 . On the technical side, arguments in 
favor of using equal weights (instead of equal heights), when applying equation 
(3) to asymmetric first-oder transitions, have been raised 29 . 

4-2 2d and 3d Ising Model 

Many real physical systems fall into the universality class of the 3d Ising model. 
Despite its simplicity, it is therefore a very rewarding object to study. Although 
many of its universal parameters have already been determined with high pre- 
cision, accurate results for some are still in the making. In particular, there has 
been recent interest in the universal surface tension and the critical-isotherm 
amplitude ratios 30 . To obtain them, one needs accurate interfacial tension 
results below the Curie temperature. Here multimagnetical (MUMA) simula- 
tions 6 have become the enabling technique for Binder's 23 histogram method, 
which was originally proposed in this context. 

For the 2d Ising model Onsager's exact result could be reproduced with 
good accuracy. However, for the 3d model the temperature dependence of 
MUMA interfacial tensions has come out fairly erratic. Therefore, the results 
of Ref. 31 should be regarded as the up-to-date best estimates. Meanwhile, 
considerable technical improvements of MUMA calculations are feasible (see 
also the next subsection) and it seems worthwhile to start off a new, large 
scale, MUMA 3d Ising model simulation. 

4.3 SU(3) Gauge Theory 

One is interested in the interfacial tension for the confinement/deconfinement 
phase transition. The use of MUCA techniques has been explored by Gross- 
mann, Laursen et al. 32 ' 33 . In particular, they noticed that it is suitable to use 
an asymmetric lattice, V — L z L 2 L t with L z > 3L. This forces the interfaces 
into the L 2 plane and ensures a flat region for the minimum of equation (3), 
thus greatly facilitating the extraction of finite-lattice values for the interfacial 
tension and, consequently, the FSS analysis. 
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For SU (3) gauge theory the interfacial tension is usually denoted by the 
symbol a and estimates are conveniently given as multiples of T 3 , where T c 
is the dcconfincment temperature. Using the conventions of , the estimates 
of 33 seem to be 

a = 0.052(4) T 3 , (L t = 2) and a = 0.020(2) T 3 , (L t = 4) 

This may be compared with the later estimate by Iwasaki et al. 34 

a = 0.02925(22) T 3 , (L t = 4) and a = 0.0218(33) T 3 , (L t = 6) 

The discrepancy (only L t — 4 can be compared) is presumably due to too small 
lattice sizes in the MUCA study. Physically, one is interested in the L t — > oo 
limit. Recently, it has been suggested that the the strong L t dependence can 
be eliminated by using tadpole improved actions. Beinlich et al. 35 report 

cr = 0.0155(16) T 3 , (L t = 2 and L t = 3) . 

4-4 Electroweak Phase Transition 

Baryon violating processes are unsuppresscd for T > T c , where T c is the elec- 
troweak critical temperature. It has been conjectured, that this may allow to 
explain the baryon asymmetry in nature. Models tie the nucleation rate to 
the interface tension of the transition. Using an effective scalar field theory 37 
or the full theory 38 , MUCA and related techniques turn out to be useful for 
simulations at a Higgs mass 

m H w (35 - 37) GeV , 

where one deals with a relatively strong first-order transitions, as needed to 
explain the baryon asymmetry. Unfortunately, it turn out that the transition 
weakens for higher Higgs masses, see Ref. 39 for a concise review. 

5 Systems with Conflicting Constraints 

In these systems one encounters large free energy barriers due to disorder 
and frustrations. MUCA simulations try to overcome the barriers through 
excursions into the disordered phase. Examples are spin glasses, proteins (see 
Hansmann in these proceedings), hard optimization problems and others. 
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5.1 Spin glasses 



MUCA studies have so far focused on the simplest, not-trivial prototypes, the 
2d and 3d Edwards- Anderson Ising (EAI) spin glass 40,41 . Significant progress 
is achieved with respect to groundstate energy and entropy calculations. How- 
ever, the slowing down with volumes size is very bad, around V 4 or, possibly, 
exponential. Certain advantages of simulated tempering are claimed in Ref. 42 . 
An uncontroversial one is that the latter approach can easily be vectorized, 
whereas other issues have straightforward (but yet untested) translations into 
the MUCA approach. In any case, presently we have no indication that simu- 
lated tempering yields a significantly improved slowing down. 

5.2 Optimization problems 

They occur in engineering, network and chip design, traffic control, and many 
other situations. General purpose algorithms for their solution are Simulated 
Annealing and Genetic Algorithms. To those, we may now like to add Multi- 
canonical Annealing 43 and Random Cost 18 . 

Multicanonical Annealing is a combination of MUCA sampling with vari- 
able upper bounds and frequent adaption of parameters. A promising study 43 
has been performed for the traveling salesman problem. Up to TV = 10, 000 
cities, randomly distributed in the unit square, were considered and scaling of 
the path length as function of N was investigated. The reported algorithmic 
performance is so good, that an independent confirmation would be highly 
desireable. 

Two non-trivial applications of the random cost algorithm emerged re- 
cently: (i) It was used for training a feed-forward multilayer perceptron, rele- 
vant for analyzing high energy physics experimental data, (ii) Topology opti- 
mization as needed for the engineering of trusses was studied 45 . 



6 Conclusions 

Sampling of broad energy distributions allows to overcome supercritical slow- 
ing down. This is well established for first-order phase transitions. Systems 
with conflicting constraints remain, despite some progress, notoriously difficult 
and for them most hope lies on achieving further algorithmic improvements. 
Finally, MUCA methods may also be of interest when dealing with second 
order phase transitions, but so far little experience exists in this direction. 
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